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Abstract 

The evolution towards equipartition in the /3-FPU chain is studied consid- 
ering as initial condition the highest frequency mode. Above an analytically 
derived energy threshold, this zone-boundary mode is shown to be modu- 
lationally unstable and to give rise to a striking localization process. The 
spontaneously created excitations have strong similarity with moving exact 
breathers solutions. But they have a finite lifetime and their dynamics is 
chaotic. These chaotic breathers are able to collect very efficiently the en- 
ergy in the chain. Therefore their size grows in time and they can transport 
a very large quantity of energy. These features can be explained analyzing 
the dynamics of perturbed exact breathers of the FPU chain. In particular, 
a close connection between the Lyapunov spectrum of the chaotic breathers 
and the Floquet spectrum of the exact ones has been found. The emergence 
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of chaotic breathers is convincingly explained by the absorption of high fre- 
quency phonons whereas a breather's metastability is for the first time iden- 
tified. The lifetime of the chaotic breather is related to the time necessary 
for the system to reach equipartition. The equipartition time turns out to be 
dependent on the system energy density e only. Moreover, such time diverges 
as £~'^ in the limit e — > and vanishes as for e — >■ oo. 



2 



Typeset using REVT^ 



I. INTRODUCTION 



In 1955, in one of the first but well known numerical simulations, Fermi, Pasta and Ulam 
(FPU) have observed the absence of thermalization in a nonlinear lattice in which the 
energy was initially fed into the lowest frequency mode. Even if a lot of progress have been 
made 0] in the last thirty years on the study of the evolution towards energy equipartition 
among linear normal modes, several points are far from being clarified. For historical reasons, 
the evolution towards equipartition has been usually analyzed considering an initial state 
where all the energy of the system was concentrated in a small packet of modes centered 
around some low frequency . Only a few studies have been devoted to the evolution from 
an initial condition where all the energy was fed in the highest frequency mode [^H. From 
these analyses it turned out that this different initial condition leads to a completely new 
dynamical behavior in the transient time preceding the final equipartition. In particular, 
the main finding was the appearance of a sharp localized mode during the transient [Q. 
Moreover, the highest frequency mode turns out to be a linearly unstable periodic solution; 
thus, starting the orbit from this initial point in phase-space put the system on a hyperbolic 
point embedded into a chaotic layer. On the other hand, when low frequency modes are 
initially excited, instability arises because of the presence of thin chaotic layers near elliptic 
points [^]. Therefore, we expect that the time-scales for the relaxation to the equipartition 
state and the physical picture of the evolution towards the final state should be strongly 
affected. 

In this paper, we present a detailed numerical and theoretical study of the dynamics of 
a /3-FPU model in the transient preceding equipartition, when only the highest frequency 
mode is initially excited. Employing the quite recent concept of breather excitations [0 , we 
are able to give a more detailed explanation of some of the behaviors observed in Ref. 0]. In 
particular, we show that the existence of such localized modes during the transient is strongly 
connected to exact breather modes for the /?-FPU model An important peculiarity of 
these excitations is that (contrary to exact breathers) they have a chaotic evolution in time, 
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therefore we have termed them "chaotic breathers" (CBs). The fact that locahzed oscillating 
excitations (that can be identified as CBs) show up spontaneously PJTO[] and persist in 



numerical simulations has suggested that they play an important role in the dynamics of 
Hamiltonian anharmonic systems. However, an important point that should be clarified is 
why these objects emerge so easily in Hamiltonian systems. Since Hamiltonian dynamics is 
reversible, large packets should break up into smaller ones at the same rate that small ones 
merge into larger one; while in the formation of a CB the latter process seems to be favored. 
It is one of the purposes of this paper to give some clarification about these points and 
in particular to emphasize the importance of these localized excitations for the transition 
towards energy equipartition. 

Another important aspect that we discuss is the existence of scaling laws for the indicators 
characterizing the approach to the equipartition state. We notice that quite general scaling 
laws indeed exist in the thermodynamic limit. Moreover, the relevant quantity for the 
equipartition time turns out to be the energy density. 

We have organized the paper in the following way. The results of numerical simulations 
concerning the appearance of CBs are presented in Sec. 0. The consequences of the existence 
of CBs for the transition to equipartition are presented in Sec. ^ together with the observed 



scaling laws. The relation between CBs and exact breathers is discussed in detail in Sec. |V 



Finally, Sec. will deal with the mechanisms of creation and destruction of breathers. Some 



final remarks and conclusions are reported in Sec. VI 



II. MODULATIONAL INSTABILITY AND ENERGY LOCALIZATION 

Denoting by M„(t) the position of the nth atom [n G [1, A^]), the equations of motion of 
the /3-FPU chain read 

Un = Un+l + Un-l - 2Un + (3 (Wn+l - - («„ - (l) 

where f3 = 0.1. The parameter (3 can be absorbed with an appropriate rescaling of Un, but 
we keep it not only for historical reasons, but also in order to make easier comparisons with 



results reported in previous papers. We have chosen periodic boundary conditions which 
allow the propagation of waves in the lattice. As recently shown in the study of the 
evolution of traveling waves with wavelength of the order of the system size (i.e. ~ A^) can 
give very useful information. In particular, the average lifetime of such waves is strictly 
related to the time necessary to reach equipartition. Here again, the propagation of waves 
and localized structures will play a fundamental role in the evolution of the system. 

As usual, in numerical simulations one does not study the real Hamiltonian system, but 
a discrete time version that approximates the time continuous dynamics. It is therefore 
essential for long time simulations to use an appropriate symplectic integration scheme in 
order to preserve as far as possible the Hamiltonian structure of the problem. We adopt 
the 6th-order Yoshida's algorithm [|12[ with a time step dt = 0.01; this choice allows us to 
obtain an energy conservation with an accuracy AE ^ 10^^^ (that corresponds to a relative 
accuracy AE/E ranging from 10^^° to 10^^^). 

We follow the approach proposed by Fermi-Pasta-Ulam in Ref. [0 where they look at the 
stability of one normal mode of the harmonic part, but contrary to them we have performed 
simulations adopting as initial condition the highest frequency mode. The highest frequency 
mode (the so called vr-mode) corresponds to the following zig-zag pattern for Un 

Un = (-1)" a and = (2) 

where a is its amplitude. Since most of the normal modes of the harmonic part of the 
Hamiltonian are no longer solutions of the full Hamiltonian, energy initially fed into one 
single mode will be shared on later times among other modes. However, this is not true for 
three particular modes that are exact solutions of the FPU lattice As the vr-mode is 
one of these solutions, in order to destabilize such initial state a small amount of noise (of 
order 10"^'*) has been added on the velocities. 

As already shown in Refs. P,p!4| the vr-mode turns out to be modulational unstable above 



a critical energy Ec, that can be analytically derived |14 



P -^^.in^t"^^ 7cos^(f)-l 
^'=-9^'^" [nJ (3cos^(f)-l)^ ■ 
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As for this mode, E = N{2a^ + 4/3a^), we easily see that if 



a > a, = sin(7r/Ar) /y^/3(9 cos2(7r/iV) - 3) (4) 

the vr-mode will be destroyed by modulational instability. 

In Fig. [l|(a) a generic evolution of the above initial state for a > Oc is reported. The grey 
scale refers to the energy residing on site n, 

En = + ^V{Un+l - Un) + {Un - M„-l) , (5) 

where the substrate potential is V{x) = + jx^. Figs. |l|(b), |l](c) and |l|(d) refer to three 
successive snapshots of the local energy En along the chain. After a very short delay, a 
slight modulation of the energy in the system appears (see Fig. |l|(b)) and the vr-mode is 
destabilized. Later, as attested by Fig. [|(a), only a few localized energy packets emerge 
from this temporary state; they correspond to oscillating localized waves and are usually 
called breathers or intrinsic localized modes. At this stage, as inelastic collisions of breathers 
have a systematic tendency to favor the growth of big breathers at the expense of small ones 
PJT5[], the number of localized objects decreases and only one very large amplitude breather- 



like excitation survives (see Fig. [I|(c)): this is the excitation we have termed CB. The CB 
moves along the lattice with a perturbed ballistic motion: sometimes the CB is even stopped 
or reflected. During its motion the CB collects energy from the visited sites on the chain, 
and its amplitude increases. It is important to note that the CB will never be at rest and 
that it propagates in general with almost the same speed (in modulus). Finally, after a very 
long time and through a mechanism we will consider later, the CB decays and the system 
reaches the equipartition of energy, as illustrated in Fig. |l](d). 

The main important aspect that arises from the above reported picture is that the 
dynamics of the model seems to favour the emergence of a well localized state (CB) during the 
transient preceding equipartition [Q. In order to give a more quantitative characterization 
of the energy localization, we introduce the following quantity 
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Co{t)=N ■ (6) 



\j=i / 

As it can be easily seen Co is of order one if Ei = E/N aX each site of the chain and of 
order N if the energy is locahzed on only one site. In Fig. ^ Cq is reported as a function of 
time for an initial condition (0) with a > a^.. Initially Cq grows in time, indicating that the 
energy, evenly distributed on the lattice at t = 0, localizes over a few sites. This localized 
state survives for some time, at later times Cq begins to decrease and finally it reaches an 
asymptotic value that is associated with a total disappearance of the CB. At this stage 
a state with a flat distribution of energy in Fourier space is attained, i.e. a state where 
equipartition of energy is fulfilled. Already at this point, it is important to remark that for 
greater than 512, all the curves Co(t) are almost coincident at any time. We discuss this 
point in Sec. 0. 

The asymptotic value Cq can be easily obtained. In the limit t — »• oo, the energy 
per site has a mean value e = E/N, but with some site dependent fluctuations, therefore 
Cq = {Ef)/{Ei)'^. Cq can be theoretically estimated within a canonical ensemble picture, 
that allows us to derive an expression for (Ef) (more details are reported in Appendix 
A). The actual Cq- value depends on energy density e and lies of course between the two 
hmiting values: corresponding to the pure harmonic case Cq = 7/4 and the pure quartic 
case Cq = 19/9 (notice that both these values are independent of e). As derived in Appendix 
A, for e = 1.44, we obtain Cq ~ 1.795, in perfect agreement with the numerical results as 
shown in Fig. ^. 

III. TRANSITION TO THE ENERGY-EQUIPARTITION STATE 

In the previous section, we have shown that the energy-equipartition state is preceded 
by a transient characterized by the emergence of a localized state. The main properties 
of this state and its relationship with exact breathers is studied in the next section. Here 



we are mainly concerned with the identification of some general aspects of the transient. 
In particular, scaling laws for two indicators measuring the relaxation to thermodynamical 
equilibrium are derived. 

As a first indicator we consider a parameter in the Fourier space S{t) that gives a 
quantitative estimation of the energy transfer among the different normal modes |TB|. The 



energy associated with the Fourier mode q = 2TTk/N, with k E {1 . . .N/2}, is in the har- 
monic approximation 0, = | (| V^p + co'(g)^lXqp), where u{q) = 2sin(g/2) represents the 
linear dispersion relation and Vq and Xq are the Fourier transforms of the velocities and 
the positions, respectively. Initially, all the energy is put in the 7r-mode, but as soon as 
the modulational instability develops also the nearest modes q = tt — 6q acquire a non-zero 
amplitude 0g, these numerical results confirm a previous analytical derivation of the most 



unstable modes [Tj]. After this initial stage, energy transfer from the highest modes to the 
lowest ones continues and the shape of the spectrum as a function of the wavevector is well 
described by an exponential 



exp [Sit) (tt - g)] (7) 



where, at any time t, S{t) is the slope of the linear fit in the Lin-Log scale. 

When the equipartition of energy is reached, the spectrum is no longer exponential (in 
other words the slope S vanishes). Therefore, the slope is a excellent tool to follow 
the transition towards equipartition. The numerical results that we have obtained yield 
the interesting conclusion that the equipartition time is a function of the energy density 
e only, as attested by Fig. H(a), where all evolutions of S(t) for various lengths are almost 



indistinguishable except for the smallest chain = 32 where finite size effects do appear [jTT 

Let us remark that in the original FPU problem, with excitations of long- wavelength 
modes, the distribution was well approximated by 0g ~ exp [— S'(t) q]. In that context. 



introducing the analytic continuation of the continuum field u{x, t) to the complex plane |]T6 
the slope was directly related to the imaginary part of the nearest singularity to the real 
axis. Work along this line would probably lead to a similar conclusion in our case. 
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The transition to equipartition can be investigated even more precisely considering the 



following indicator |18 



2 -J2 Piit) Inpg(t) 
Vit) = j^e 1 , (8) 

where Pq(t) is the probability to have an energy (pg associated to the mode q at time t; i.e. 

P.it) = ^ . (9) 
<? 

Indeed, r]{0) = whereas max(?7) = 1 is reached ||19| in the energy equipartition state where 



all = 2/N (for a detailed derivation of the equipartition value f] see Appendix B). Thus 
energy sharing among normal modes will be detected by an increase of i]{t), which can be 
considered as the percentage of modes with significant energy. It is exactly what we obtain 
in Fig. ^(b) for various chain lengths. The main conclusion that we can draw from this result 
is that again the equipartition time is a function of energy density e only, and therefore it is 
finite in the thermodynamic limit, although it may diverge in the limit of vanishing energy 
density (see the following). 

Simulations for various energy densities 0.33 < e < 10^ allow us to study the variation 
of the equipartition time as a function of the the energy density. In the low energy limit 
(0.33 < e < 1.44), a good data collapse is obtained if the time is rescaled by a typical time 
scale r ~ e"^ (see Fig. ^ (a)). At high energy (e > 1000) the typical time scale is r ~ 
as it is clearly shown in Fig. |^ (b). Let us also stress that e was also found to be the relevant 



dependence of the largest Lyapunov exponent in the equipartition state pO|Jl^ . Moreover, 
the maximal Lyapunov exponent shows two scaling laws at high and low energy density that 
coincide with those here reported for the equipartition time (namely, for its inverse) pT| , p!4| . 
However, while for the Lyapunov the transition from one regime to the other is observed at 
e ^ [pT] , p!^ , for the equipartition time the transition occurs at higher energy : e ~ 100 



for j3 = 0.1. Power law divergences of the relaxation time have been recently reported for 
several models of nonlinear oscillators and for different classes of initial conditions 
although the large limit has not been studied as carefully as in the present paper. 
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IV. RELATIONSHIP BETWEEN EXACT AND CHAOTIC BREATHERS 



In the present section, we would like to pay attention to the localized objects observed 
during the transient state, and study their relationship with the already known localized 
solutions of nonlinear lattices, that we call exact breathers. We have seen that the energy 
contained in a CB (and thus its frequency) can be very high. For example, in chains of 
N = 512 sites, we found breathers with frequency uj ~ 3.5 (when the maximal frequency 
associated to the phonon band is 2). We would like to emphasize here that the usually 
adopted Nonlinear Schrodinger approximation p5|,^ is not appropriate to describe these 
objects. 

The unstable vr-mode gives rise spontaneously to a localized breather-like excitation. 



qualitatively very similar to the exact breathers obtained for the FPU chain |j26| . The most 
evident difference is that the former move in an erratic way in the lattice and have a finite 
lifetime, while the latter are exactly periodic and mainly static. Therefore it is natural to 
try to determine to what extent it is possible to compare them. 

Another important difference between the self-created excitation and the exact breathers 
concerns their dynamical stability, which has been investigated with the aid of a Lyapunov 
analysis. On one hand exact breathers are linearly stable and remain perfectly unchanged 
during very long simulations. On the other hand, we show here that CBs are strongly 
chaotic. 

In order to compute the Lyapunov exponents of the system, we have used the standard 
algorithm of Benettin et al . However, we are here interested in the Lyapunov exponents 
at short times instead of the asymptotic ones. We evaluate the cumulative average corre- 
sponding to such exponents as a function of time. In Fig. ^ the evolution of the localization 
parameter Co(t) is reported (panel (a)) together with the corresponding cumulative average 
for the first four Lyapunov exponents (panel (b)) for a typical evolution starting from an 
unstable vr-mode. 

It is clear from Fig. |^(b) that the localized object, which has spontaneously appeared, 
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is chaotic because its presence (indicated by the peak in Co{t)) is associated with a positive 
maximal Lyapunov exponent. Moreover, another typical feature of the CB is that the 
running average for the maximal Lyapunov exponent has an higher value during the transient 
than in the equipartition state. The evolution of this running average is quite similar to that 
of Co(t): there is an initial growth followed by a decrease at later times. Naively, one would 
usually expect that an higher degree of chaoticity should be related to a higher degree of 
energy equipartition: the present simulations clearly shows the opposite. 

One should stress that the running average for the Lyapunov exponents relaxes very 
slowly to the asymptotic value, because its value is affected by the estimation at earlier 
times. A more efficient method to obtain the asymptotic Lyapunov exponents is therefore 
to restart the running averages after the transition to equipartition state (i.e. after the 
decrease of Cq). The dash-triple dotted line of Fig. ^(b) shows that the convergence is much 
faster if this is done. 

Another interesting result is the fact that a gap is present in the distribution of the Lya- 
punov exponents during the transient: the first Lyapunov is clearly above the others whereas, 
when the system reaches the energy equipartition state, the gap disappears (as shown in 
Fig. ^(b)) and the whole spectrum of Lyapunov exponents is approximately linear |^ . 



Furthermore, a measure of the localization of the Lyapunov vectors shows that during 
the transient the first Lyapunov vector is localized contrary to all others. This localization 
of the Lyapunov vector disappears after that the equipartition is reached, since the energy 



density is below the strong stochasticity threshold, see . An interesting question is to see 
whether this peculiar structure of the tangent space (i.e. the gap in the spectrum and the 
localization of the first Lyapunov vector) can be related to the properties of exact breathers. 

Let us briefiy recall that exact discrete breathers are time-periodic and spatially expo- 
nentially localized solutions of the equations of motion (|I]). Their frequency is always higher 
than the frequency of the top of the phonon band. Calculating such a solution is equivalent 
to looking for a fixed point of the stroboscopic T map 
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{{Un, Mn})(0) ^ T{{Un, Un}) = {Un, Un}(tb), (10) 

where th is the period of the solution. This can be done with the aid of a Newton process, 
using the tangent map dT. The latter relates linearly an initial perturbation {ej, ej}(0) to its 
image {ei,ei}{tb) where the perturbations evolve according to the N linearized equations 
of motion 

'in = [l + 3/?(m„+i - M„)^] (e„+i - e„) - [l + 3/?(m„ - (e„ - e„_i). (11) 

Starting from a sufficiently good approximation, the Newton method converges to a periodic 
solution, satisfying Eq. ( p!0D up to machine precision pGf. One can then investigate the 
linear stability of the breather solution with a standard Floquet analysis, i.e. computing 
the eigenvalues of the 2N x 2N matrix dT (a periodic solution is linearly stable when all 
eigenvalues lie on the unit circle of the complex plane). 

The main results concerning exact breathers in the FPU chain can be summarized as 
follows. They exist for every frequency above the phonon band; moreover the spatially 
antisymmetric solutions (centered between two particles, as the generic example plotted in 
Fig. ^(a) and sometimes called P-modes [^) are linearly stable, while the symmetric (the 
S'T-modes BO]) are unstable. The spectrum of the Floquet matrix, discussed in detail in 



31| , consists of a "continuum" of spatially extended eigenmodes (the linear phonons) and 
a discrete part with a spatially symmetric and exponentially localized mode. Due to the 
time reversibility of the solution, the real axis is a symmetry axis of the spectrum (if A is 
an eigenvalue of dT, then A* is also an eigenvalue). Fig. ^(b) shows schematically half of 
the spectrum (the other half is its complex conjugate) with the continuum and the discrete 
localized mode. As the problem is solved for a finite chain, the continuum corresponds in 
reality to (A^ — 1) modes with a higher density of modes close to the end of the band. 

Such linear modes out of the phonon band and with opposite symmetry to the original 
solution, have been observed and studied in the framework of Klein-Gordon chains. They are 
called pinning or translation modes . Indeed, an excitation of the periodic solution 



in the direction of this mode produces an oscillation of the center of energy or can even lead 
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to a propagation of the solution along the chain without radiative decay. We have checked 
that these features are also present in the /5-FPU model. 

The problem now is to link these exact solutions with the self-created excitations which 
have been observed in the simulations. The essential difference is that the former are lin- 
early stable while the latter are strongly chaotic. However the comparison of the Lyapunov 
spectrum of the CB and the Floquet spectrum of the exact breather solution shows strong 
similarities: they both present a "continuum" of extended states and an isolated spatially 
localized mode. A quite natural hypothesis is that non-zero perturbations of the exact 
breather could lead to a strong destabilization of the translation mode. If this is true, 
one can understand why this instability is not destructive for the breather: a perturbation 
along this direction only leads to a coherent displacement of the breather and in addition, 
it explains why, in the simulations, the non-linear excitation is very seldom static. 

In order to verify that the most unstable direction corresponds to the translation mode of 
the exact solution, we have performed the following test. We study the evolution of an exact 
solution, once it has been perturbed with a spatially antisymmetric gaussian noise all over 
the lattice. The antisymmetry of the initial condition allows the breather to remain at rest 
centered between two sites and this simplifies the comparison. The Lyapunov analysis reveals 
that the perturbed state is now chaotic. The maximal Lyapunov eigenvector associated with 
the perturbed solution and the Floquet vector corresponding to the translational mode for 
the unperturbed case are shown in Fig. |^. The two vectors are both symmetric and in 
perfect agreement. This is a convincing evidence that the structure of the tangent space of 
a CB can be qualitatively interpreted in terms of the Floquet spectrum of "neighboring" 
exact periodic solutions. 

We would like to stress that the fact that breather-like excitations move is not a sufficient 
condition for their chaoticity. Similarly to what has already been found numerically for 
Klein-Gordon chains [^], many exact mobile breather solutions can be exhibited by the 
FPU model (see pi|-p^] or for a discussion on approximations of moving breathers). 
Such solutions satisfy the relation 
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{Un,Un}iT) = {Un-l,Un-l}{0) (12) 

where T, the inverse of the velocity, is a multiple of the period of internal vibration of the 
breather. One should also note that a small amplitude phonon tail dresses these solutions, 
such that one can consider that the breather is in equilibrium with the emitted and the 
absorbed radiation (for more details see Ref. [Q). Many of these solutions are linearly 
stable and were checked to have zero maximal Lyapunov exponent; the chaotic regime is 
reached only if the whole system is sufficiently perturbed. The direction of the translation 
mode is easily excitable, but is not the direction responsible for chaoticity: thus mobile 
breathers are not necessarily chaotic and strong chaos occurs only when other modes are 
sufficiently excited. 

V. FORMATION AND DESTRUCTION OF BREATHERS 

At this stage, we have understood the connection of the CB with exact FPU breathers 
but we should explain how the destabilization of the vr-mode can lead to only one localized 
solution. Indeed, as attested by Fig. |l|(a), the modulational instability gives rise first to a 
few packets of energy but then, because of their interactions, energy is concentrated in only 
one. 

A. The localization process 

In order to study this effect in a controlled manner, we have put four exact breathers on 
the lattice: three moving ones with small amplitude and a frequency u = 2.12, that is just 
above the phonons' upper band edge (u; = 2), and one at rest with a larger amplitude and a 
frequency uj = 2.75. As shown in Fig. due to collisions, the biggest breather successively 
absorbs the smaller ones and gives rise to a moving large amplitude breather-like excitation. 
Such a result is a generic example of the collision process in the FPU chain and it explains 
why a single localized CB emerges during the transient. 
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Let us emphasize that even if some locahzation processes have aheady been reported 
in homogeneous nonhnear lattices 0, the process presented here is particularly interesting 
because of the absence of the self-regulation process. Indeed, in the Klein-Gordon systems 
studied in Ref. P,p!0|, the process is regulated by a stronger pinning effect. Discreteness 



provides a path to localization but is also responsible for the pinning effect. This is stronger 
for larger excitations (big breathers are easier trapped), and this does not allow a collapse 
of all the breathers into a single very large excitation. For the FPU model, pinning effects 
are perhaps too small to be detected. Moreover, we have observed that the velocity of the 
CB slightly increases with its energy. As a consequence, we see in Fig. |^ the emergence of 
only one large amplitude breather-like excitation from the initial vr-mode. 

Is this effect present also in the thermodynamic limit ? As the lifetime r and the velocity 
f of a CB are finite, one could predict that, for very long chains, the CB will probably not 
have enough time to collect all the energy present in the system. If, in small chains, the 
breather could easily move through the system many times, in very long chains this would 
not be the case. One can expect that above a critical chain length L^. = vt more than 
one CB will be observed in the transient preceding equipartition. To check this, we have 
investigated the case with an initial amplitude a = 0.8 for the vr-mode; in this situation the 
typical velocity of the final breather is f ~ 0.2 and its lifetime is r ~ 10^. From these data 
the critical length Lc is ~ 2 ■ 10^. For chains > 10^ and a = 0.8, additional simulations 
have shown that this is indeed true: the energy is no longer localized in only one huge CB 
but in few of them. 

A long chain could therefore be considered as a juxtaposition of almost non-interacting 
sub-chains of length Lc. And in each sub-chain, one single CB is created by the mod- 
ulational instability. This interpretation is consistent with the observed saturation of 
the localization parameter Cq when the size of the system increases (Fig. 0). Indeed, 
a chain can be considered as made up of 2 independent parts A and B, and one has 
Co{A + B) = {Cq{A) + Co{B))/2 = Co{A) (Co is an intensive-like quantity). As discussed in 
Sec. m, the lifetime is a function of the energy density e. Thus the critical length of course 
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depends on e. 

One can define two different regimes associated to the growth of the excitation. In a first 
stage the main breather present in the chain absorbs the smaller ones by colhsions. After this 
initial stage the growth slows down. This effect could be explained by an observation made 
by Bang and Peyrard for the Klein-Gordon equation ||T5|]. They noticed that the energy 
transfer between localized modes depends on their relative energy difference: the transfer is 
less efficient when the energy difference (or frequency difference) increases. For our system 
we have observed that, after the first collisions, the lattice contains one big breather with 
a rather high frequency and a large population of high frequency phonons waves. The 
presence of the latter is due to the destabilization of the original 7r-mode via modulational 
instability Correspondingly, the energy transfer-rate between the phonons and the CB 
is still positive but reduced with respect to the initial one. We want to stress that this second 
regime was detected because of the very small initial perturbations of the vr-mode, since for 
bigger perturbation this regime can easily be missed |Q. In addition, in this paper we have 
concentrated our attention on the energy region just above the critical energy defined by 
Eq. (^, where the FPU chain evolves towards equipartition on a very long time scale. 

B. Breather's metastability 

Once the localization process is finished, a very large amplitude breather is moving in 
the system and all the high frequency waves have been absorbed. The hungry excitation 
cannot grow any more and this regime corresponds to the plateau where Co has reached its 
maximum. At later times, the excitation disappears and now we want to explain how this 
phenomenon can take place. We have seen that the effect of collisions of the CB with high 
frequency modes leads to absorption of the latter. During such a scattering process, small 
quantities of energy radiate towards low frequency modes. Therefore, after the initial stage 
the only phonons left in the system are the low frequency ones. It is natural to suspect that 
the destruction of the CB is related to their interaction with the low frequency modes. 
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In order to understand this interaction, we have performed numerical experiments where 
one exact breather initially at rest collides with the following wave packets of low frequency 
phonons centered on the site Uq: 

Unit) = A cos(gn — ut) e ^'o . (13) 

We have taken q = 0.2, i = 13 and an amplitude A G [0.5; 3]. Fig. ^ shows that the 
collision process is very efficient in destroying the CB excitation. 

As shown in Fig. ^(a), the decrease of the energy of the breather Eb can be considered as 
linear at the earlier times with a good approximation. Moreover, a more detailed study of the 
phenomenon shows that the amplitude of the wave packet determines the destruction rate. 
In particular, the slope of the decrease of Eb is an exponential function of the amplitude. 
This result suggests that a very complicated nonlinear mechanism is at the origin of the 
interaction. To give a more quantitative explanation we should study the scattering on one 
phonon with the CB as done by Cretegny et al ||31]] for exact breathers at rest in the Klein- 
Gordon model. However, here, the study will be much more complicated because small 
perturbations will easily put the FPU-breather in motion. 

A possible explanation of the CBs metastability is the following: at the beginning of the 
simulation, when only the vr-mode is excited, and even during the modulational instability 
process, the low frequency phonons are not present. During the growing process of the CB 
the low frequency phonon band is populated due to radiative processes. When all the high 
frequency phonons have been absorbed by the CB, the excitation can only loose energy 
due to the destructive interaction with the low frequency phonons. Therefore, we can say 
that the CB "die of starvation" . The destruction of the CB is associated with a significant 
increase in the population of low frequency linear waves. This transition corresponds to 
the final decrease of the slope S{t) plotted on Fig. ^ and to the final relaxation to the 
equipartition state. 
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VI. CONCLUSION 



We have seen that the evolution towards equipartition in the /3-FPU chain, starting 
from the 7r-mode as initial condition, gives rise to a striking localization process. The 
spontaneously created excitations are moving breather-like excitations with a finite lifetime 
and a chaotic dynamics. The features of these transient localized modes can be explained 
by exploiting the correspondence with exact breathers. 

We want also to stress that, contrary to the previous belief, we have found that the 
localization mechanism originally obtained in the Klein-Gordon case |^,n]|, is also a transient 
state: the breathers obtained by collisions have a very long lifetime but are indeed CBs as 
characterized by the local Lyapunov exponent. The lifetime of the process is however much 
longer in this case, because pinning effects on the excitation are stronger. 

It is important to recall that the FPU equation was at the origin of the rediscovery of the 



soliton in the continuum limit, which justified the recurrence phenomenon observed by 
Fermi-Pasta- Ulam [|l| (see also [^). In the present paper, we have shown the fundamental 
role played by a second family of excitations, the breathers, in achieving the equipartition 
state on a discrete FPU lattice. The creation of a CB can be considered as an efficient 
mechanism to transfer energy from high frequency modes to low frequency phonons. 

Let us also remark that the breathers are created from the high frequency zone corre- 
sponding to zero group velocity. This portion of the phonon band was also recently found 



to be crucial for the transition to equipartition ||3^. Indeed, the transition to equipartition 
disappears if a renormalized FPU Hamiltonian with a dispersion relation without zero group 
velocity region is considered |3^. As in electromagnetism and string theory, the idea that 
this ingredient is absolutely necessary to observe the transition to equipartition of energy is 
another reason to justify a posteriori our study. 

Moreover, we have found that the fundamental parameter for the dynamics of the system 
is always the energy density e = E/N (i.e. energy /volume). This is true not only for the 



value of the maximal Lyapunov exponent (as already noticed in Refs. 0,50]) but also for the 
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transient time towards energy equipartition. One should of course emphasize that this is the 
only parameter which makes sense in the thermodynamic limit. Therefore, the transition to 
equipartition happens in a finite time in the thermodynamic limit, but this time diverges as 
an inverse power of the energy density in the zero energy density (zero temperature) limit. 

As a final remark, this study emphasizes that the concept of a breather is not only 
important for the energy localization in lattices, but is also crucial to address one of the 
main classical problem of statistical mechanics: the transition to equipartition. 
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FIGURES 

FIG. 1. Evolution of the local energy En along the chain. In Fig. (a), the horizontal axis 
indicates the position along the chain and the vertical axis corresponds to time (time is going 
upward). The grey scale goes from En = (white) to the maximum E'n-value (black). The lower 
rectangle corresponds to < t < 3000 and the upper one to 5.994 10^ < t < 6.10^ Figs, (b), (c) 
and (d) show the instantaneous En along the = 128 chain at three different times. Note the 
difference in vertical amplitude. The initial vr-mode amplitude was a = 0.4 (while Oc = 0.0317). 

FIG. 2. Evolution of Co{t) for chains of various lengths: N = 32, 128, 512, 1024, 2048 and 
8192. Each curve corresponds to the average over 20 simulations with the same energy density 
e = E/N = 1.44 but with different random initial noise added to the velocities. For all the reported 



simulations E > Ec 111]- 



FIG. 3. Fig. (a) (resp. (b)) presents the evolution of the slope S{t) (resp. r]{t)) versus time 
for chains of various lengths; =32, 64, 128, 256 and 512. Each curve corresponds to the average 
over 20 different initial conditions with the same energy density e = 1.44. 

FIG. 4. In Fig. (a) (resp. (b)) the evolution of ri{t) versus the rescaled time t (^^^ (resp. 
t [j^j ) is reported for different energy densities. Each curve in Fig. (a) (resp. (b)) corresponds 
to the average over 20 (resp. 50) different initial conditions for a chain with N = 512 sites. The 
dot-dashed line reported in Fig. b indicates the asymptotic value f] ~ 0.655, which is theoretically 
estimated in Appendix B. 

FIG. 5. Fig. (a) presents the evolution of Co(i) and Fig. (b) presents the cumulative average 
of the first four Lyapunov exponents when the initial condition is a vr-mode with a = 0.4 in a 
A^ = 128 chain. The solid line corresponds to the first Lyapunov exponent (the dash-triple dotted 
line corresponds to the value of the largest Lyapunov exponent after the transient), the dotted one 
to the second, the dashed one to the third and the dash-dotted line to the fourth. t2 and in (a) 
are defined in Fig. 1. 
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FIG. 6. Fig. (a) presents the spatial shape of an exact hnearly stable FPU breather, solution 
of Eq. (1) with a frequency uj = 2.5. Fig. (b) shows the schematic repartition of its Floquet 
eigenvalues on the complex plane. 

FIG. 7. Shape of the translation mode. The solid line corresponds to the eigenmode of the 
Floquet analysis whereas the stars correspond to the first Lyapunov eigenvector. We have reported 
only the part of the vectors relative to the positions when the internal phase of the breather is zero 
(i.e. when the kinetic energy is minimal). 

FIG. 8. Breathers' merging. Fig. (a) shows the energy evolution of an exact breather initially 
at rest {u = 2.75) after collisions with three small exact moving breathers {iv = 2.12). Fig. (b) 
(resp. (c)) shows the energy repartition at t = (resp. at t = 2 • 10"^). 

FIG. 9. Breathers' destruction. Fig. (a) shows the evolution of the energy of an initial exact 
static breather (a; = 2.4) after collisions with a wave packet of phonons {A = 1.5, Iq = 13, q = 0.2). 
Fig. (b) (resp. (c)) shows the initial (resp. final) repartition of energy. 



APPENDIX A: CALCULATION OF THE EQUIPARTITION VALUE Cq OF THE 

LOCALIZATION PARAMETER 



Introducing the usual parameter (3 = I/Ze^T, we obtain directly from Eq. 
totic expression corresponding to energy equipartition 



Co = iV. 



{Ed' 



where the spatial averages are 



{E^) 



K Via 
'2 2 2 



the asymp- 



(Al) 



(A2) 



and 
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For long chains, we can evaluate {Ei) and {E'f) using the canonical ensemble. Introducing 
the configurational partition function 

/ + CO 
exp [-(5V{x)] dx , (A4) 
-oo 

we have (F,) = -if and {V^) = 
One obtains therefore: 

and 

In the pure harmonic case this gives (Ei) — ^ and (Ef) — whereas in the pure 
quartic case we obtain (Ei) = ^ and (Ef) = 

In conclusion, = | in the pure harmonic case and Cq — ^ ~ 2.11 in the pure 
quartic case. But one can even compute the result using the complete FPU-potential. 

V{x) = ^ + S^. 

Using fCi the modified Bessel function of the second kind, we have 

4 

We obtain finally for 5 = 0.1, the value 

Co ~ 1.795 , (A8) 
in excellent agreement with the numerical value. 
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APPENDIX B: ESTIMATION OF THE EQUIPARTITION VALUE fj IN THE 

HARMONIC APPROXIMATION 



We would like to give a theoretical estimation of the indicator ri{t) (see Eq. (^)) in 
the limit t — »• oo, i.e. in the equipartition state (the asymptotic value is indicated as 
lim ri{t) = f]). 

t — *oo 

Let us rewrite the spectral entropy as 

N/2 

S = -Y. P.W lnp,(t) (Bl) 

9=1 

where Pq{t) = . Introducing the energy per mode e = — ^ 0g, in the equipartition 

09 ^ 9=1 

9 

state, we can easily derive the following expression 

5- = Z — Z ^ine + lnfyj (B2) 

2 

The average appearing in (|B2D is < 01n0 >= -T7^^4>q^^4>q^ due to the equipartition 

9 

the index q has been neglected. 

The corresponding expression for the quantity (|) is 

r^ = e e--.«t>'^'t» (B3) 

Assuming that the fluctuations can be neglected, we obtain 

< 01n0 >=< >< ln0 >=< > In < >= elne (B4) 

and fj would be exactly one. However, the numerical value of f/ estimated in the low energy 
limit is quite different from one (namely, it is ~ 0.795). This indicates that the inclusion of 
the fluctuations is fundamental to give a realistic estimate of f/ (this was firstly noticed in 
). In order to consider the fiuctuations, let us write 

0" - 1 

ln0 = lim . (B5) 

n^O n 



This allows one to reexpress Eq. (p^) as a function of average of powers of 
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< > -e 

< In >= lim —- (B6) 

n^O n 

The estimation of fj is now reduced to the estimation of terms of the type < 0" >, which 
can be again evaluated within the canonical ensemble, after a mode inverse temperature f3* 
is introduced. In the harmonic approximation 

/ IT/ 12 _|_ , ,21 |2\ " n „| 



with 



A=/"dJ^)'e-^-#=r(*+l)^^, (B8) 



oo 



where F is the Euler Gamma function. In particular, e =< >= (/?*) ^ and the expression 
( P?!) reduces to 

< 0" >= e"r(n + 1) = e"n! . (B9) 
We are now able to give an expression for the spectral entropy. 

< 01n0 >= lim ^ ^ ' — - = e[log(e) + r'(2)] , (BIO) 

being r'(2) = 1 — 7, where 7 ~ 0.5772 is the Euler constant. Finally we obtain 

fj = e^-i ~ 0.655 . (BU) 



This results have been also obtained with a different approach in [0. We notice that this 



value is in perfect agreement with the numerical one found in the high energy limit (see Fig. 
4(b)), while it is an underestimation of the one found in the low energy limit (see Fig. 4(a)). 

The origin of such a discrepancy is related to the structure of the Hamiltonian when 
expressed in terms of normal modes. Due to the quartic term in the potential, the mode 



interaction matrix is no longer diagonal |T3|. However, in the high energy limit a sort 



of "random phase approximation" should be valid for the terms appearing in the mode 
interaction matrix and the off-diagonal terms should average to zero. Therefore we are left, 
in the action-angle representation {Iq,6g), with diagonal terms of the form 
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H,, = {u;{q)Q + ^uj{qfll < cos' 9, > . (B12) 

In the limit N ^ oo the nonhnear term in Hgq becomes neghgible and the harmonic ap- 
proximation is recovered. 

In the low energy limit, the action-angle variables are strongly correlated and a "ran- 
dom phase approximation" is no longer valid. Therefore the off-diagonal terms cannot be 
neglected and give a contribution (although comparatively small) to fj. 
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FIGURE 3 Cretegny, Dauxois, Ruffo, Torcini 
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